Obtaining models of subsurface properties is of value to extractive industries (for example, hydrocarbon production or mining) Direct sampling via drill cores can provide constraints, but direct sampling is infeasible over large areas, particularly in challenging environments. Remote sensing via geophysical data (e.g. seismic waves, gravity, electromagnetic waves) is used instead to develop property models over large areas. However, inversions of remotely sensed geophysical data are typically non-unique, that is, a range of geologic models is consistent with the measured data. Jointly inverting independent geophysical data sets has been proposed as a method to reduce ambiguity in the resulting property model (e.g., Vozoff and Jupp, 1975).
The concept of joint inversion is well known in the geophysics community.
Generally, a penalty function (sometimes called an objective function, a cost function, or a misfit function) is formulated as a function Fj of observed geophysical data dj and predicted data dj*(m) from the candidate subsurface model for each data type:Ψ=λ1F1(d1,d1*(m)+λ2F2(d2,d2*(m)+ . . . +R(m)  (1)
The goal of the inversion is to determine the model m that minimizes equation (1). R(m) is a regularization term that places constraints on the model, for example, model smoothness or length. The variables λi represent arbitrary weights that determine the relative contributions of the data types to the penalty function. In general, the inversion is non-linear and it is difficult to obtain a solution that satisfies equation (1) directly. Common solutions to this problem are to perform iterative inversions, using gradient based approaches or stochastic approaches (e.g., Monte Carlo or genetic algorithms).
A typical joint inversion algorithm is disclosed by Moorkamp et al. (2010), and is shown schematically in FIG. 1. Geophysical data 11 are collected over a region of interest. An initial property model 12 is estimated, and simulation data 14 is generated from the model. A joint penalty function 13 is constructed along the lines of equation (1) based on the measured data and simulated data. Convergence is tested at step 15 by determining whether the penalty function indicates an update is possible. If an update is possible, the property model is updated at step 16, and the inversion loop begins a new iteration. If an update is no longer possible, the inversion is said to have converged to a final property model 17 (or failed if the data is not sufficiently fit).
Geophysical data (FIG. 2) may include gravity data, electromagnetic data (magnetotelluric or controlled source), seismic data (active or natural source), or other types of data. Property models (FIG. 3) may include one or more rock properties, such as fluid saturation, porosity, resistivity, seismic velocity, density, lithology, or other properties.
The vast majority of joint inversion algorithms use constant weight values for the duration of the inversion (e.g. Lines et al., 1988, Johnson et al., 1996, Julia et al., 2000, Linde et al., 2006). Unfortunately, as weight choices are arbitrary, this may direct the solver through model space in an inefficient manner, resulting in an increased number of iterations to achieve convergence. Further, in many problems model space is populated with many models that are local minima of equation (1). These local minima can “trap” gradient based solvers, leading to incomplete convergence. This is diagramed in FIG. 4. Finally, if constant weights are used, several different inversion runs must be completed using different weights in order to find an optimal solution.
Colombo et al. (US patent application publication No. 2008/0059075) disclose a method for joint inversion of geophysical data and applications for exploration. However, as described above, their method is prone to trapping in local minima Tonellot et al. (US patent application publication No. 2010/0004870) also disclose a method for joint inversion of geophysical data, but their method uses constant weights on data terms as well. Both of these techniques would benefit from a new method to mitigate local minima.
Lovatini et al. (PCT patent application publication No. WO 2009/126566) disclose a method for joint inversion of geophysical data. The authors use a probabilistic inversion algorithm to search model space. These types of solvers are not hindered by local minima, but are much more computationally intensive than the gradient-based solvers to which the present invention pertains.
Publications in a variety of different fields describe methods to adaptively change weights. For example, in U.S. Pat. No. 7,895,562, Gray et al. describe an adaptive weighting scheme for layout optimization, in which the importance of a priority is scaled based on the magnitude of a lesser priority. Unfortunately, this method does not allow the priority weights to be adjusted during the solution of the layout optimization problem.
In U.S. Pat. No. 7,487,133, Kropaczek et al. describe a method for adaptively determining weight factors for objective functions. In this method, the data component with the highest penalty (i.e., the data component that contributes most to the objective function) receives an increased emphasis (weight) in a subsequent penalty function. However, simply increasing emphasis on the component with highest penalty can result in misleading results if that component is trapped in a local minimum. In this situation, increasing the weight of this component will ultimately result in a converged result. Though other components are not satisfied, they are down-weighted in the penalty function and therefore no longer contribute to the final model.
Chandler presents a joint inversion method in U.S. Pat. No. 7,383,128 using generalized composite weight factors that are computed during each solver iteration. These weight factors are related to the independent variables (i.e., data). In this method the weights are chosen in such a way as to render error deviations to be represented by a non-skewed homogeneous uncertainty distribution. However, this method requires the input of a priori error estimates and is limited to two-dimensional problems, and is not applicable to joint inversion of multiple independent data sets because the derivation of this method is limited to variables that can be represented as “orthogonal coordinate-oriented data-point projections,” which is not the case for geophysical data sets.
Adaptive weights have also been considered in applications to neural network optimization algorithms. For example, Yoshihara (U.S. Pat. No. 5,253,327) discloses a method in which synaptic weights are changed in response to the degree of similarity between input data and current synaptic weight. This emphasizes connections with high importance in the network, resulting in a more efficient search of model space. Unfortunately this and similar disclosures (for example, Jin et al., U.S. Pat. No. 7,363,280, or Ehsani et al., U.S. Pat. No. 5,917,942) actually are “adaptive learning” systems, in which past searches of model space are used to guide the future search of model space. While this can help avoid local minima and can result in more efficient computation in some situations, it depends on historical experience in the solution space, rather than intrinsic properties of the data themselves.